# compute the mean of a bunch of numbers...
from pyblog import *
import math

@var_dist
def mean(): return Normal(0, 1e-6)      # weak prior

@var_dist
def prec(): return Gamma(1, 1)

@var_dist
def coord(i): return Normal(mean(), prec())

def generate(numpoints):
  return query([mean(), prec()] + [coord(i) for i in range(numpoints)], [],
               burnin=0, scans=0, stats=False, outtyp = QUERY_LAST)

def inference(data, scans):
  ans = query([mean(), prec()], [coord(i) == d for i,d in enumerate(data)],
              burnin = int(.1 * scans), scans=scans, stats=True)

  return [x.mean() for x in ans]

if __name__ == "__main__":
  data = generate(100)
  print "True mean, std. dev:", data[0], math.pow(data[1], -.5)
  ans = inference(data[2:], 1000)
  print "Inferred mean, std. dev:", ans[0], math.pow(ans[1], -.5)

  ## # for fun see what happens without the lifted variables optimization
  ## configure("PARAM_LIFTED_VARIABLES", False)
  ## ans = inference(data[2:], 1000)
  ## print "Inferred mean, std. dev:", ans[0], math.pow(ans[1], -.5)
  ## 
  ## # re-enable the optimization it should still work
  ## configure("PARAM_LIFTED_VARIABLES", True)
  ## ans = inference(data[2:], 1000)
  ## print "Inferred mean, std. dev:", ans[0], math.pow(ans[1], -.5)
  
  
